Dynamical simulation of current fluctuations in a dissipative two-state system 
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Current fluctuations in a dissipative two-state system have been studied using a novel quantum dynamics 
simulation method. After a transformation of the path integrals, the tunneling dynamics is computed by deter- 
ministic integration over the real-time paths under the influence of colored noise. The nature of the transition 
from coherent to incoherent dynamics at low temperatures is re-examined. 

PACS numbers: 02.70.Lq, 05.30.-d, 05.40.+j 



A two-state system coupled to a dissipative environment 
is the archetypical model for tunneling phenomena in con- 
densed phase. It has found widespread applications in solid- 
state physics [jl|-[|], most recently in interlayer charge trans- 
port in high-r c superconductors [Q| and quantum computing 
|gP, as well as in biophysics [||] for the modeling of electron 
transport in biochemical reactions. One of the most intrigu- 
ing features of this model is a dynamical phase transition be- 
tween coherent tunneling and incoherent relaxation. This was 
first predicted by Chakravarty and Leggett [^||] and later con- 
firmed by experiments on interstitial tunneling in niobium |^]. 

Although the existence of the coherent-incoherent transi- 
tion is widely accepted, its precise nature and location has 
been called into question by some recent calculations [10-13|]. 
Coherence is a phenomenon of dynamics, yet an exact treat- 
ment of tunneling in the time domain has so far been out of 
reach. The original prediction of the transition [j|] was based 
on a dynamical but approximate theory, whereas the more re- 
cent theories, suggesting the transition would occur at a much 
weaker damping than predicted earlier, were based on statisti- 
cal mechanical calculations [ 10,0], 



In this Letter, we describe a new exact numerical method 
for calculating the real-time dynamics of dissipative quantum 
systems and use it to investigate the transition from coher- 
ent to incoherent dynamics in a two-state system at low tem- 
peratures. Previously, the only exact numerical approach to 
tunneling dynamics has been the dynamical quantum Monte 
Carlo (QMC) method UQ. But all real-time QMC simu- 
lations fail at longer times because the signal-to-noise ratio 
of the results vanishes exponentially due to the highly os- 
cillatory integrand. This problem is commonly referred to 
as the dynamical sign problem. In the case of large band- 
width of the dissipative environment, QMC simulations fur- 
ther suffer from a slowing-down problem caused by the in- 
creasingly long-lived correlations in the sampling process. 
The new method eliminates both problems through a gener- 
alized Hubbard-Stratonovich transformation and allows us to 
perform the functional integration over paths of the tunneling 
system by a deterministic method, while statistically sampling 
fluctuations from an ensemble of Gaussian noise trajectories. 

Dissipative two-state systems are often described by the 
spin-boson model, 



+q^C v {a v + al), 



(1) 



where q = q^a z /2 is the position operator of the tunnel- 
ing system with intrinsic tunneling frequency A, a x and a z 
are Pauli spin matrices, and h = 1. The effect of the har- 
monic environment is fully characterized by a spectral den- 
sity J(u>) = tt Yl u C*3 — u k )> f° r which the Ohmic form 
J(u>) = 2-Kauj/q^ is experimentally the most relevant and 
theoretically the most interesting. The Ohmic spectral density 
introduces a single dimensionless damping constant a. This 
model must be regularized by an upper cutoff uj c of the spec- 
tral density. The scaling limit ui c ^> A is characteristic of 
tunneling in solids and as shown by scaling arguments [jOfl, 
the Ohmic spin-boson model has nontrivial dynamics only for 
a < 1, and the renormalized tunneling frequency 



A r = A (A/w c ) 



a/(l- 



(2) 



is the only frequency scale of the dynamics at zero tempera- 
ture other than uj c . The transition from coherent to incoher- 
ent dynamics occurs at a critical damping a c at which the Q 
factor of the tunneling oscillations vanishes. Tunneling oscil- 
lations can be observed as a damped oscillatory component 
of the position correlation function C(t) = | (cr z (t)a z (0) + 
&z(0)o~ z (t)) that is present in addition to an incoherent re- 
laxation background Jlrj], The asymptotic long-time behav- 
ior is always dominated by an algebraic incoherent decay, 
C(t) oc aA-H- 2 0. 

A coherence criterion equivalent to finite Q is a finite de- 
phasing time of the quantum beats that manifest themselves 
as tunneling oscillations. A measure of this dephasing time is 
given by the lifetime of delocalized states of the tunneling sys- 
tem |J[] (see also |+) ± which are eigenstates 
of the tunneling current j — Aa y /2 — & z /2. The correlation 
time r of the current correlation function 



c jj (t) = o(j'(*)j(o)+i(o)j(<)> 



(3) 



equals the lifetime of these superposition states. A finite 
correlation time of the current correlation function thus im- 
plies coherent oscillations in the position correlation function. 
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This relationship allows us to identify coherence even for very 
strongly damped cases in which oscillations may be masked 
by the incoherent background. For a two-state system, the re- 
lation C(t) = —ACjjit) provides another direct connection 
to previous studies on the position correlation function. 

To compute the exact dynamics of Cjj, we employ a 
hybrid stochastic/deterministic numerical method which we 
shall label chromostochastic quantum dynamics (CSQD). 
This method, which is generally applicable to quantum sys- 
tems with linear dissipation, is based on the path integral for- 
mulation of dissipative quantum dynamics The time 
evolution of the reduced density matrix for the system coor- 
dinate q can formally be represented by a double functional 
integral ^ 

fit rq't 

p(q l ,q' i ;t) = / V[q] / V[q']e iSo[q] - iSo{q] F[<l^]- (4) 

So[q] is the action of the undamped quantum system, and 
its interactions with the environment are incorporated into a 
complex-valued influence functional F[q,q'] — exp(— $' — 
with 

<&W]= fdt' fdt"{q{t')-q'{t')) 

Jt Jt 

xL'(t'-t")(q(t")-q'(t")), (5) 

&'[q,q'} = \ fdt' fdt"{q{t')-q'(t')) 

xm 1 (t'-t")(q(t")+q'(t")). (6) 

L'{t) is the real part of the autocorrelation function of the col- 
lective bath mode J2 U C u (a v + a\ ), averaged over an ensem- 
ble of free oscillators, m is the mass of the tunneling parti- 
cle, and j(t) is the classical friction kernel associated with 
the spectral density J(ui). This influence functional describes 
a 'factorized' initial preparation, i.e. with the particle con- 
strained to one side for times t < to with the environment 
fully relaxed. To simulate equilibrium correlation functions, 
it is necessary to push the preparation back to a sufficiently 
large negative time to < [ pOQ and insert measurement oper- 
ators j at times and t. 

A direct quantum Monte Carlo evaluation of a time- 
discretized version of type (Q) can be prohibitively expensive 
due to the dynamical sign and slowing-down problems. An 
alternative approach for single-particle problems or dissipa- 
tive systems near the Markovian limit is the explicit iteration 
of a short-time propagator. This reduces the path integral to a 
manageable series of matrix multiplications J22|]. Another al- 
gorithm was recently presented by Cao, Unger and Voth [ ^3| ] 
for an environment with a few oscillators. By directly sam- 
pling the oscillator paths and propagating the system coordi- 
nate for each sample, they eliminate memory effects. 

In this Letter, we tackle the problem of memory effects in 
the important case of a continuum of environmental modes, 
such as in the case of Ohmic friction, where sampling over 



individual oscillator trajectories is not feasible. The major ob- 
stacle in trying to decompose the path integral (Q) into short- 
time propagators lies in the interaction kernel L'(t) because 
its range diverges as the temperature approaches zero. In com- 
parison, the friction kernel ^(t) poses no problem because it 
vanishes at time larger than uj~ 1 , the shortest timescale in our 
problem. 

The problem with the long-range interactions introduced by 
L'{t) can nonetheless be solved, albeit at the cost of introduc- 
ing an additional path variable. The exponential of the non- 
local action $'[q, q 1 ] can be decomposed into a superposition 
of time-local phase factors, 

exp(-$W]) = 

J mw[t]ew{-ij*?mw)-<fV))Y (7) 

The distribution W[£] is real and Gaussian, with 
(€(t)£(t'))w = L'(t - t'), and normalizable through the 
condition $'[(7,5] = 0. Formally, this decomposition is a 
Hubbard-Stratonovich transformation in a function space over 
the interval [to , t] . Equation (^) is equivalent to the construc- 
tion of an influence functional for a classical colored noise 
source [fTT]], and as such, we will interpret the function as 
a noise trajectory. 

In the CSQD algorithm, the propagation of the system co- 
ordinates q and q' is carried out deterministically for each re- 
alization of the noise trajectory, while the noise trajectories 
are sampled from the distribution W[£\. Instead of generat- 
ing weights from a Metropolis-type random walk, statistical 
weights are assigned to noise trajectories a priori by numeri- 
cally filtering white noise. 

In general, there are two ways to treat the remaining term 
If q acts on a finite-dimensional Hilbert space, q{t) con- 
tains discrete jumps between the eigenvalues of q. The 'state 
vector' that is propagated must then remember these 'virtual 
transitions' during a finite number of preceding time slices. 
But since the friction kernel ^(t) decays rapidly over the 
memory time w" 1 , the number of transitions needed 'in mem- 
ory' is small for large cutoff. On the other hand, if q describes 
a continuous degree of freedom, the kinetic term in So [q] re- 
quires the relevant paths q(t) to be smooth. Then the limit 
uj c — > 00 can be applied to (g|), making $" [q, q'] a time-local 
functional. 

For the two-state system, we achieve excellent accuracy 
with a maximum of two tunneling events 'in memory'. The 
error resulting from this truncation can be made arbitrarily 
small by increasing uj c , i.e., moving further into the scaling 
regime p4| ] . The length of our time steps is a fixed fraction 
of the inverse bandwidth uj~ 1 , which is short enough to make 
the particular choice of approximate short-time propagator a 
minor issue. 

Empirically, we find that the statistical variance of CSQD 
grows no faster than logarithmically with t. This remark- 
able performance compared to the exponential divergence ob- 
served in conventional QMC methods results from the fact 
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that deterministic integration is not affected by the oscillatory 
nature of the integrand. Consecutive samples in the CSQD 
simulation are by design statistically independent. This com- 
pletely eliminates any slowing-down problem. We have tested 
the validity and accuracy of our method, comparing to known 
analytic results for the real-time fluctuations and the relax- 
ation behavior of the tunneling coordinate at damping a — 
1/2 and a <C 1. We find excellent agreement between our 
numerical results and previous analytic calculations. 

Fig. [l] shows CSQD results for Cjj(t) at zero temperature 
and for < a < 1/2. Except for the initial drop at times 
t < w" 1 , a ll f° ur curves represent functions of the scaling 
argument A r t only. To obtain a formal estimate for r, we 
define it to be the first zero of Cjj(t) J2~5|]. Positive current- 
current correlations do persist up to times of order A^ 1 , i.e., 
the scaled dephasing time A r r is nonzero and finite (at least) 
up to a = 0.4, although it declines significantly with increas- 
ing a. Quantum coherence becomes increasingly short-lived, 
however, making it impossible to resolve oscillations from a 
background of monotonous signal decay for a > 0.3 (Fig. ||). 

An expanded view of the initial decay of Cjj (t) in the vicin- 
ity of the critical value a c = 1/2 is provided in Fig. ||. In 
this region, the correlation time r decreases rapidly with in- 
creasing a (note the different scales on the time axes). Com- 
paring results for different cutoff frequencies uo c , we observe 
markedly different scaling behavior for a above or below 1/2. 
For a < 1/2, r is scaled by A r , leaving A r r finite; whereas 
for a > 1/2, r scales as i.e. A r r vanishes in the scaling 
limit. 

Critical behavior at a c = 1 /2 is also predicted by perturba- 
tion theory in the bare tunneling frequency A. In the scaling 
limit, the short-time behavior of the current correlation func- 
tion is given by 

A 2 

Cjj(t) = -± cos(ira)(A r t)- 2a . (8) 

For a < 1/2, this expression (valid for u^ 1 C t < A) as 
well as the exact perturbative result (valid for t A -1 ) is 
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FIG. 1. Current correlation function for various damping 
a < 1/2, Lu c /A r = 50. 



0.5 


y i i i 

1 OC = 0.4 " 


0.0 








0.0 
0.0 


\ 




-0.5 






I I I I I I I I I I L. 

2.0 4.0 6.0 ^ 



FIG. 2. Coherent oscillations and long-time decay of the current 
correlation function (u! c /A r = 50). 

positive. It follows that r > A -1 , which can only be satisfied 
if t scales with A r rather than u> c . For a > 1/2, (||) is nega- 
tive, and we conclude that Cjj (t) turns negative at a time r of 
the order of uj~ 1 . 

Taken together, these results lead to a clear picture of 
how the scaled dephasing time A r r depends on the damp- 
ing strength, shown in Fig. |j. We see a gradual decline of 
A r r to very small but finite values, followed by a discontin- 
uous drop at a c = 1/2. The finite value A r r ps 0.034 at 
a = 1/2 emerges both from our simulations and from ap- 
plying Cjj(t) = —C(t)/A to analytic results for the position 
correlation function at a = 1/2 [ |16[|26| |. 

Our conclusion a c = 1/2 coincides exactly with the value 
originally found by Chakravarty and Leggett [Q] for the dy- 
namics of a two-state system with a factorized initial prepara- 
tion. In contrast, Lesage, Saleur, and Skorik JIo| ] as well as 
Costi and Kiefer [|ll]] recently found a lower value a c « 1/3 
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FIG. 3. Scaling behavior of the current correlation function in the 
vicinity of its first zero. 
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for the critical damping. In both their works, the criterion of 
coherence was not the observation of tunneling oscillations 
themselves, but the presence of inelastic peaks in the spectral 
function of the position fluctuations. For strongly damped os- 
cillations in the symmetric function C(t), however, the inelas- 
tic peaks at positive and negative frequencies become wide 
enough to merge into a single peak centered at uj = 0. The 
absence of a peak at finite frequency therefore does not neces- 
sarily indicate the absence of coherence. In order to determine 
the true critical coupling, a more elaborate analysis needs to 
be performed on the spectra of [ |To|Jl~T] ] . 

In conclusion, we have introduced a new numerical algo- 
rithm for the dynamics of dissipative quantum systems, which 
solves the dynamical sign problem and eliminates slowing- 
down problems. Its validity, efficiency and accuracy have 
been demonstrated for the spin-boson model with Ohmic dis- 
sipation. We have determined the lifetime of coherent su- 
perposition states at zero temperature and different strengths 
of the damping from the current correlation function. This 
timescale decreases with increasing damping strength a, but 
remains finite (indicating quantum coherence) for all a below 
the critical value a c — 1/2. Analytic results, where available, 
are reproduced with good agreement. Details of the CSQD 
method and its applications to more complex problems such 
as the noise spectrum of fractional quantum Hall systems will 
be reported elsewhere. 
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